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Abstract 

We construct doubly periodic Stokes flows in two dimensions using 
elliptic functions. This method has advantages when the doubly periodic 
lattice of obstacles has less than maximal symmetry. We find the mean 
flow through an arbitrary lattice in response to a pressure gradient in an 
arbitrary direction, and show in a typical example that the shorter of the 
two period lattice vectors is an "easy direction" for the flow, an eigenvector 
of the conductance tensor corresponding to maximal conductance. 

It is known, and we rederive it below, that two-dimensional (2D) Stokes flows 
can be represented in terms of two complex analytic functions [U It is plausible, 
then, that doubly periodic 2D Stokes flows should have a representation in 
terms of doubly periodic complex analytic functions, that is, elliptic functions 
0E1- Such a representation was promised in 1959 by H. Hasimoto 0], but 
it did not appear. Other authors in the succeeding decades alluded to such a 
representation [3], and even, like Hasimoto, quoted results following from it 
The Hasimoto article may have appeared much later as lecture notes in Japanese 
\7\. Meanwhile, there are other, perhaps more straightforward, ways to represent 
doubly periodic 2D Stokes flows. These include matching of flows around a 
single obstacle across periodic cell boundaries |5J E3 EH , integral equation 
methods |12j . biharmonic solvers on a grid and finite element methods. 

In spite of this long history, we have thought it useful to present the elliptic 
function method, because there is still, apparently, no readily available descrip- 
tion of it. Furthermore, as we shall show, this approach solves one aspect of the 
problem which is not at all simple in the most common cell matching approach, 
namely the appropriate boundary condition for flow through a general periodic 
lattice in a general direction. With this method we describe the typical flow 
through a generic lattice. 
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1 2D Stokes Flows 

Let us represent the 2D flow with velocity vector field u(x, y) as a complex scalar 
function u by means of the usual isomorphism 



u = u T x 



u v y <-> u(x,y) = u x + iu y (1) 



Regarding the x-y plane as the complex z plane, we note that the divergence 
and curl of u are given by 

divu = 2$t(du/dz) = (2) 
curl u = 2^(du/dz) = uj (3) 

Here the first equation expresses incompressibility of the flow u, and the second 
defines the vorticity lo, understood as the (scalar) component of a vector field 
normal to the plane. Similarly, the gradient of the pressure P in this complex 
representation is 

VP <-» 2dP/dz (4) 

The Stokes equation, representing the balance of viscous stress in the fluid by 
the pressure gradient, then becomes 

/iV 2 u = VP <-> 4^id 2 u/9zdz =2dP/dz (5) 

where /i is the viscosity. In light of Eqs. and © this means that 

d{P/ti-iu})/dz = (6) 

that is, that 

f(z) = P/fi - ilo (7) 

is a single- valued function of z, holomorphic except for possible poles 1 . If we 
choose such a function /, we can integrate Eqs. J2J and © to find the flow 

u = ~(zf-J fdz + g) (8) 

where g{z) is a second holomorphic function, having logarithmic singularities at 
the poles of /. Thus u is represented in terms of two holomorphic functions, / 
and g. 

The force exerted by the flow u on a finite obstacle, given by closed contour 
C, just involves the enclosed residues of /. The force on a small, oriented line 
segment Az due to the fluid on its right is 

AF = iPAz + 2ifi(du/dz)Az (9) 

Now, using Eqs. (JZJ) and (jHJ), and integrating over the closed curve C, oriented 
in the conventional positive direction, we find the force on C due to the fluid 
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outside it 



la 

= f 

y 



(f + f)dz+f {zf' + g')dz 
c Jc 



d{zf + g) + f fdz 
c Jc 



ifx <j> f dz — —2ir[i res(/) 



The last line follows because u in Eq. © is single-valued. 



(10) 

(11) 
(12) 



2 Pressure in Doubly Periodic Flows 

Consider the integer lattice generated by two complex numbers, ui\ and lo 3 , with 
^s(uj 3 /lui) > 0, consisting of the points 



W mn = 2mwi + 2nUJ 3 



(13) 



for all integers m and n. Suppose identical obstacles are located at these places, 
forming a doubly periodic array. A Stokes flow through this array, represented 
as in Eq. JSJ, is characterized by a function / with very restrictive properties. 
It is single- valued, it is holomorphic outside the obstacles, it has poles inside 
the obstacles, its imaginary part is doubly periodic, and its real part is a doubly 
periodic function plus a linear function, where the linear function is essentially 
the average pressure < P >, increasing linearly along its (constant) gradient. If 
we further ask for the simplest function of this type, having only simple poles 
at the W mn , then there is essentially only one possibility, the Weierstrass zeta 
function £(z), with a linear correction. This follows from the theory of elliptic 
functions 0- (Note that C depends also on the lattice constants lo\ and 013, but 
we regard these as fixed parameters and do not indicate this dependence.) 
More precisely we make use of the quasi-periodicity of C( z )j 



C(z + 2cu a ) = C(z) + 2 Va 
for a = 1,3, where rj a — C,{uj a ) are constants satisfying 

771 u> 3 - r/sujx = iir/2 

Then we can take, for the function /, either of 

— 1 

fa = 1 A^aC(z) ~ T) a z\ , 

where a — 1,3. Using Eqs. (fH}> and l(TB|) . we verify that 

A(« + 2wi) - h(z) 

f 1 (z + 2uj 3 ) = fi{z) - 7r/|wi| 

h{z + 2ui) = / 3 (*)+7T/|w3| 

h(z + 2u 3 ) = f 3 (z) 



(14) 



(15) 



(16) 



(17) 
(18) 
(19) 
(20) 



3 



Thus in both functions only the real part, which we interpret as P/n, shows a 
linear growth. The average pressure gradient corresponding to f a is perpendic- 
ular to tu a . From its projection on the other lattice vector we determine that it 
is 

<VP> <- 2 < dP/cfz >= - . (21) 

2|w Q |3(u;iW3) 

By taking linear combinations of fx and we can find a suitable / with the 
pressure gradient in any direction with respect to the lattice. The combination 
corresponding to a pressure gradient of the same magnitude as fx and f% but in 
the direction e %B is 

/ = e i6 ( " Cz (22) 

with 

'7?i3?(e- ifl w 3 ) - ^(e^V 



C = i U V ^ _ V " (23) 
V 3(wi^ 3 ) / 

The force on the fluid in a period parallelogram D, with sides 2u>x and 2^3, 
due to the average pressure in the flow corresponding to f a , is 

F P = f <P>idz = i < dP/dz > dzdz (24) 

JdD J JD 

= -2<dP/dz> [ [ dxdy^^^ (25) 
J JD \^a\ 

The force on the fluid due to each obstacle in this same flow, since C,(z) is 
normalized to have residue 1 at each simple pole, is —2iL0 a im/\u) a \, by Eq. (|12|) . 
There is on average one obstacle in each period parallelogram. Thus the force on 
the fluid due to the average pressure (applied somehow from outside) is balanced 
by the force due to the obstacles, and the net force on the fluid is zero, as is 
always the case in Stokes flows. 

A more general function / in the representation of Eq. JSJl for the Stokes 
flow u would be a superposition of translates of fx and fa, always keeping 
the singularities inside the obstacles, and not in the physical region of flow. 
Equivalently, one could take a multipole expansion of such functions. This 
would be a series in derivatives of fx and f^,. These functions are doubly periodic 
with higher order poles at the lattice points, that is, they are the Weierstrass V 
function and its derivatives. Thus / has the form 

N c 

f(z)^e te C-Cz + Y,c n r in) (26) 

n 

If the boundary contours have reflection symmetry in the origin, then it is 
enough to take only the terms in the series with n odd. 
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3 Velocity Field in Doubly Periodic Flows 



From Eq. (jSJ) we know that the flow u corresponding to / above is 

1 / _ _ W = \ 

u = - ze- % \ - C\z\ 2 - 2e w In \a\ - Cz 2 /2 + z £ c„7>(«) _ ^ c n V [n ~ X) + 5 

\ n 11 / 

(27) 

where g is a second analytic function, still to be determined. We have already 
used the freedom to add terms of the form g in adding a term e l6 lner, where a is 
the sigma function, another of the special functions of elliptic function theory [2] . 
Such logarithmic terms compensate the multivaluedness of logarithmic terms in 
the integral of /. 

Now among all the functions g that we could choose, we want the one that 
makes u doubly periodic, and that makes u — (say) on the boundary contours 
of the obstacles. Since the logarithmic part of g has been explicitly written as 
a separate term, g as defined above is single- valued, and hence has a Laurent 
series in a neighborhood of the origin. We must anticipate, though, that g has 
singularities at all the lattice points W mn , and hence that this series would con- 
verge only out to the nearest such lattice point. This difficulty can be removed 
because double periodicity of u determines the principal part of g at all lattice 
points. We see, for example, that u blows up logarithmically at z = but, be- 
cause of the poles in it blows up like 1/z at every other lattice point. This is 
inconsistent with double periodicity, and hence g must cancel the 1/z behavior 
at all W mn except 0. Also, in the Laurent series for g about z = 0, there will be 
negative powers of z. Translates of these terms contribute to the principal part 
of g at other W mn . Let S be a set of indices (m, n) sufficiently large to label all 
W mn in a disk large enough to contain a period parallelogram. Then we take g 
in the form 

W N » N " 

9=- £ j±w~+Y, h i £ (s-w^-'+xv ( 28 ) 

t \ ran - n , . „ ... 

(m,n)GS 3=0 (rn,n)eS 3 = 1 

Because the poles of g have been included as explicit terms, the power series 
with coefficients aj now converges out to the nearest W mn not indexed in 5*, 
and we choose this index set large enough to get good numerical behavior in 
our final step. 

Lastly we specify the boundary conditions that determine the coefficients a n , 
b n , and c„ in Eqs. I|27|) and l|28|l . We choose a region D with a boundary dD 
that is entirely within the fluid region, and such that every point z in dD has 
a corresponding point z' in dD displaced by one of the lattice constants ±2cj q . 
Typically D would be the period parallelogram centered on the origin with 
vertices (±^,±^3). We ask that the complex function u be periodic, taking 
the same value at points z and z' in dD related by a lattice constant, and also 
that the derivative du/dz be periodic in the same sense. We also require that u 
take prescribed values on the boundary contour C of the obstacle. Choosing a 
large number of pairs of points (z, z') on dD and also points z on C for imposing 
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these conditions, we obtain an overdetermined inhomogcncous system of linear 
equations that we then solve in the sense of least squares. Because c n occurs 
together with its complex conjugate c^ t the system must be regarded as real 
linear, not complex linear. With suitable choices for the truncation parameters 
N a , Nb, and N c , this boundary value problem has a solution accurate to many 
decimal places. We confirm the solutions in [S], for flows through square and 
triangular lattices of obstacles, to 5 decimal places, typically. 

It is superficially surprising that one can stipulate conditions on both u 
and its derivative in this elliptic problem, but it is clear physically that such 
periodic solutions must exist. It is also clear mathematically, if one considers 
the equivalent Stokes flow on a torus. This argument only applies, though, if the 
applied stress that drives the flow has all the properties that we have required of 
the function / in Section [2J For other functions / the boundary value problem 
would have no solution. Knowing the form of / in advance, as we do, the actual 
boundary values of P and u emerge as part of the solution. Other approaches 
to this problem are typically restricted to situations in which one knows the 
boundary values by symmetry, where P is constant on part of dD 7 for example, 
or where u is normal or tangential on 3D, special cases that do not hold for 
general lattices and flows. 

Let us define the mean of a flow u to be the constant flow V such that 
the flux of V and the flux of u into a period parallelogram are the same. By 
incompressibility of the flow, this is also the flux of V (and u) out of a period 
parallelogram. If we take the period parallelogram to be the one with vertices 
(±cji, iws), then we have the formula 



V 



23(^1^) 



3 ( I udz ] u>i + I / udz ) u>3 



(29) 



The integrals depend only on the endpoints, by incompressibility of the flow, 
and the formula is invariant under a common translation of the endpoints. Thus 
the choice of period parallelogram is in fact arbitrary. 

Using this formulation we have computed the tensorial relation between the 
force F in Eq. (|12l) and the mean velocity V in Eq. 129fl in many examples. For 
regular square and triangular lattices, F and V are simply proportional, but 
more generally 

V = GF (30) 

where the conductance tensor G is a real symmetric matrix (and we are now 
representing V and F as real 2- vectors, not complex scalars). A typical com- 
puted flow is shown in Figure [T] Here io\ = 1, 0)3 = 0.2 + 0.5i, and the circular 
obstacles have radius 0.2. The pressure gradient F is in the y direction, but 
the mean flow V is approximately along 0)3, the shorter of the lattice vectors. 
The reason for this behavior is simple, and in accord with common sense. The 
direction of 0J3 is an "easy direction" for the flow, an eigenvector of G, and 
the corresponding conductance eigenvalue for the rather wide "channels" along 
this direction is large. The obstacles, being relatively closely spaced along the 
direction of W3, create "walls" along the sides of the channels. The conductance 
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Figure 1 : A Stokes flow through a lattice of circular obstacles is shown in four 
period parallelograms. The lattice vectors are u>\ = 1, u>s = 0.2 + 0.5i, and the 
obstacles have radius 0.2. The pressure gradient is vertical, but the mean flow 
is, to good approximation, along w 3 , the shorter of the lattice vectors. 

corresponding to the other, perpendicular eigenvector, normal to the "walls", 
is much smaller than the conductance along the channels. The ratio of the two 
eigenvalues, a measure of the anisotropy of the conductance tensor, is about 
2.545 in this example, noticeably more anisotropic than the lattice itself. It is 
clear that the conductance anisotropy could be indefinitely large for this lattice 
if the round obstacles were larger and almost touched along the W3 direction. 

The theory described here provides a starting point for problems involving 
flow near ciliated surfaces. Within a ciliated layer, the cilia, hairlike projections 
from the surface, would constitute the regular array of obstacles. In the simplest 
case they would be rigid and perpendicular to the surface. Consider, for exam- 
ple, a straight pipe of circular cross section, with such a ciliated inner surface, 
still leaving a free cylindrical channel down the center. Suppose the length of 
the cilia is much greater than the mean distance between them, and much less 
than the radius of the pipe. The Poiseuille flow problem asks for the flow in 
response to a pressure difference between one end of the pipe and the other. 
The usual boundary condition says that the flow should be zero on the walls of 
the pipe, and this condition still determines the flow in a thin region near the 
wall, comparable in thickness to the mean spacing of the cilia, but within most 
of the ciliated layer, the mean flow will be given in the manner described here. 
Via a similarly thin transition layer, the Poiseuille flow down the unobstructed 
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central channel matches this mean flow at its edge. Systematic anisotropy in 
the placement of the cilia, including, one must anticipate, a merely statistical 
anisotropy, would impart rotation to the fluid channel, as it amounts to a "ri- 
fling" of the inner surface of the pipe. Considerations of this kind are relevant 
in biological settings where such surfaces are very common. The cilia as we 
have just described them are merely passive obstacles, but often they are active 
agents driving flows. The methods of this paper could also be a starting point 
for treating such phenomena, as we hope to show in the future. 
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